Pocillopora Juvenile Fusion
Ariana Huffmyer (R Markdown author)
2019
Setup
Set up workspace, set options, and load required packages.
rm(list=ls(all=TRUE)) knitr::opts_chunk$set(root.dir = "~/R",warning=FALSE, message=FALSE)library(gridExtra)
library(multcomp)
library(emmeans)
library(ggplot2)
library(dplyr)
library(effects)
library(plyr)
library(lattice)
library(glmmTMB)
library(effects)
library(lsmeans)
library(MuMIn)
library(MASS)
library(car)
library(FSA)
library(lmerTest)
library(lme4)
library(blmeco)
library(ggsci)
library(coxme)
library(survival)
library(cowplot)
library(RColorBrewer)
library(tidyverse)
library(partitionBEFsp)
library(factoextra)
library(ggfortify)
library(cowplot)
library(survminer)1. GROW OUT SURVIVORSHIP
Load data
Survivorship is measured as binary on each individual spat in the experiment.
#load data
rearing<-read.csv("Data/GrowOut_Responses_PostGenotyping.csv", header=TRUE, sep=",", na.strings="NA")
rearing$Quadrant<-as.factor(rearing$Quadrant)
rearing$Diversity.Type<-relevel(rearing$Diversity.Type, ref="Individual")
#rearing$Total.Juveniles<-as.factor(rearing$Total.Juveniles)
#rearing$Fusion.Partners<-as.factor(rearing$Fusion.Partners)
#rearing$Richness<-as.factor(rearing$Richness)Subset the final timepoint for additional analyses.
final<-rearing[which(rearing$Timepoint == "Final"),]
multiple<-rearing[which(rearing$Community == "Multiple"),]
multiple_final<-multiple[which(multiple$Timepoint == "Final"),]Analysis
Analyze with a logistic binomial regression model.
Analyze at the final timepoint.
#with final dataset
survmod2<-glmer(Survivorship ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Colony) + (1|Tank) + (1|Slide), data=final, family=binomial)
summary(survmod2)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula:
## Survivorship ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners +
## (1 | Parent.Site/Colony) + (1 | Tank) + (1 | Slide)
## Data: final
##
## AIC BIC logLik deviance df.resid
## 224.9 255.5 -103.5 206.9 212
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -8.6175 -0.6293 0.3049 0.5167 1.4715
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 3.186e-01 0.5644532
## Colony:Parent.Site (Intercept) 1.488e-02 0.1219734
## Parent.Site (Intercept) 1.634e-07 0.0004042
## Tank (Intercept) 1.709e-01 0.4133998
## Number of obs: 221, groups:
## Slide, 115; Colony:Parent.Site, 21; Parent.Site, 6; Tank, 5
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.450767 0.519055 2.795 0.00519 **
## Richness -0.007283 0.267024 -0.027 0.97824
## CommunityMultiple -1.847623 0.618356 -2.988 0.00281 **
## Fusion.Partners 3.036311 0.946317 3.209 0.00133 **
## Richness:Fusion.Partners -0.432505 0.319527 -1.354 0.17587
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Rchnss CmmntM Fsn.Pr
## Richness -0.503
## CmmntyMltpl -0.231 -0.556
## Fusn.Prtnrs -0.105 0.414 -0.505
## Rchnss:Fs.P 0.145 -0.474 0.405 -0.902
## convergence code: 0
## Model failed to converge with max|grad| = 0.00237872 (tol = 0.001, component 1)Anova(survmod2)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Survivorship
## Chisq Df Pr(>Chisq)
## Richness 0.5759 1 0.447925
## Community 8.9279 1 0.002808 **
## Fusion.Partners 21.1261 1 4.3e-06 ***
## Richness:Fusion.Partners 1.8322 1 0.175871
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1#final dataset
summary(effect(c("Richness", "Fusion.Partners"), xlevels=4, survmod2))##
## Richness*Fusion.Partners effect
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.5201784 0.9361096 0.9949754 0.9996265
## 2 0.5183603 0.9041990 0.9880627 0.9986242
## 3 0.5165418 0.8587551 0.9719082 0.9949464
## 4 0.5147229 0.7966036 0.9353234 0.9816172
##
## Lower 95 Percent Confidence Limits
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.3287517 0.8160149 0.9478899 0.9855863
## 2 0.3643003 0.7954607 0.9392540 0.9827277
## 3 0.3149904 0.7068390 0.8811552 0.9520535
## 4 0.2323719 0.5215122 0.6410654 0.7004323
##
## Upper 95 Percent Confidence Limits
## Fusion.Partners
## Richness 0 1 2 3
## 1 0.7058580 0.9797581 0.9995363 0.9999905
## 2 0.6690068 0.9581692 0.9977482 0.9998920
## 3 0.7128515 0.9387678 0.9938440 0.9994880
## 4 0.7879786 0.9336587 0.9915323 0.9991807summary(effect("Fusion.Partners", xlevels=4, survmod2))##
## Fusion.Partners effect
## Fusion.Partners
## 0 1 2 3
## 0.5186977 0.9110372 0.9898281 0.9989197
##
## Lower 95 Percent Confidence Limits
## Fusion.Partners
## 0 1 2 3
## 0.3642350 0.8027598 0.9426089 0.9840155
##
## Upper 95 Percent Confidence Limits
## Fusion.Partners
## 0 1 2 3
## 0.6696673 0.9626407 0.9982685 0.9999280summary(effect("Community", xlevels=4, survmod2))##
## Community effect
## Community
## Individual Multiple
## 0.9480743 0.7421160
##
## Lower 95 Percent Confidence Limits
## Community
## Individual Multiple
## 0.8373723 0.5918728
##
## Upper 95 Percent Confidence Limits
## Community
## Individual Multiple
## 0.9847894 0.8509763dispersion_glmer(survmod2) #no overdispersion## [1] 0.9357489plot(residuals(survmod2)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)Create plot for effect of community type.
eff.surv1 <- predictorEffect("Community", survmod2)
surv.effects1<-plot(eff.surv1,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bars", width=4, col="black"),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="right",
xlab=expression(bold("Number of Juveniles")),
main=""); surv.effects1
#lattice=list(key.args=list(space="right",
#border=FALSE,
#title=expression(bold("Genotypic \nRichness")),
#cex=1,
#cex.title=1)));surv.effects1Create plots for fusion partner effects.
survmod2a<-glmer(Survivorship ~ Richness + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Colony) + (1|Tank) + (1|Slide), data=final, family=binomial, subset=c(Community=="Multiple"))
eff.surv2 <- predictorEffect("Fusion.Partners", survmod2a, xlevels=4, rug=TRUE)
surv.effects2<-plot(eff.surv2,
lines=list(multiline=TRUE,
col=c("lightblue","mediumblue", "blue4", "darkblue"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="right",
xlab=expression(bold("Fusion Partners")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));surv.effects2
survmod2b<-glmer(Survivorship ~ Fusion.Partners + (1|Parent.Site/Colony) + (1|Tank) + (1|Slide), data=final, family=binomial, subset=c(Community=="Multiple"))
eff.surv2b <- predictorEffect("Fusion.Partners", survmod2b, xlevels=4, rug=TRUE)
surv.effects2b<-plot(eff.surv2b,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="right",
xlab=expression(bold("Fusion Partners")),
main="",
lattice=list(key.args=list(space="none",
border=FALSE,
title=expression(bold("")),
cex=1,
cex.title=1)));surv.effects2b
Generate mean values of survivorship for 1) number of juveniles and 2) within “multiple” juvenile groups, mean values for number of fusion partners.
meansurv <- ddply(final, c("Community"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meansurv## Community N mean sd se max lower upper
## 1 Individual 58 0.7758621 0.4206554 0.05523477 1 0.3552066 1.196518
## 2 Multiple 163 0.6687117 0.4721270 0.03697984 1 0.1965847 1.140839meansurv2 <- ddply(multiple_final, c("Fusion.Partners"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meansurv2## Fusion.Partners N mean sd se max lower upper
## 1 0 75 0.3466667 0.4791133 0.05532324 1 -0.1324466 0.825780
## 2 1 48 0.9583333 0.2019409 0.02914766 1 0.7563924 1.160274
## 3 2 24 0.9166667 0.2823299 0.05763034 1 0.6343368 1.198997
## 4 3 16 0.9375000 0.2500000 0.06250000 1 0.6875000 1.1875002. GROW OUT FUSION
Analysis
Analyze with a logistic binomial regression model.
Final timepoint dataset.
Analyze fusion as a product of genotypic richness within “multiple” juvenile groups.
#with final dataset
fusemod1<-glmer(Fusion ~ Richness + (1|Parent.Site/Colony) + (1|Tank) + (1|Slide), data=final, family=binomial, subset=c(Community=="Multiple"))
summary(fusemod1)
Anova(fusemod1, type="II") # no effects on rate of fusion by the end point
dispersion_glmer(fusemod1) #no overdispersion
plot(residuals(fusemod1)) + abline(h=0, lty=2)
There is an effect of genotypic diversity on fusion rates. Plot the relationship.
eff.fuse1 <- predictorEffect("Richness", fusemod1, xlevels=4, rug=TRUE)
fuse.effects1<-plot(eff.fuse1,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Prob(Successful Fusion)")),
xlab=expression(bold("Richness")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));fuse.effects1
Generate a summary of proportion fused at the end of the grow-out period.
meanfusion <- ddply(multiple_final, c("Richness"), summarise,
N = length(Fusion[!is.na(Fusion)]),
mean = mean(Fusion, na.rm=TRUE),
sd = sd(Fusion, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Fusion, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanfusion## Richness N mean sd se max lower upper
## 1 1 59 0.5423729 0.5024778 0.06541704 1 0.03989509 1.0448507
## 2 2 52 0.7115385 0.4574670 0.06343925 1 0.25407150 1.1690054
## 3 3 28 0.1428571 0.3563483 0.06734350 1 -0.21349118 0.4992055
## 4 4 24 0.6250000 0.4945354 0.10094661 1 0.13046464 1.11953543. GROW OUT GROWTH
Analysis
Growth is measured as polyp expansion in each spat by number of spat in each sibling type and number of fused individuals. Growth is measured in final dataset as all spat started with 1 polyp each.
Analyze polyp expansion with the final dataset and present means.
polypmod1<-glmer(Polyps ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners + (1|Parent.Site/Colony) + (1|Tank) + (1|Slide), data=final, family = poisson)
summary(polypmod1)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: poisson ( log )
## Formula:
## Polyps ~ Richness + Community + Fusion.Partners + Richness:Fusion.Partners +
## (1 | Parent.Site/Colony) + (1 | Tank) + (1 | Slide)
## Data: final
##
## AIC BIC logLik deviance df.resid
## 590.6 618.1 -286.3 572.6 148
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.1240 -0.5469 -0.1021 0.4563 1.4273
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 1.858e-09 4.311e-05
## Colony:Parent.Site (Intercept) 7.237e-02 2.690e-01
## Parent.Site (Intercept) 1.445e-03 3.802e-02
## Tank (Intercept) 2.789e-03 5.282e-02
## Number of obs: 157, groups:
## Slide, 92; Colony:Parent.Site, 21; Parent.Site, 6; Tank, 5
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.279359 0.147887 8.651 <2e-16 ***
## Richness -0.021274 0.091952 -0.231 0.817
## CommunityMultiple -0.095335 0.160455 -0.594 0.552
## Fusion.Partners 0.068683 0.143007 0.480 0.631
## Richness:Fusion.Partners -0.009974 0.052539 -0.190 0.849
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Rchnss CmmntM Fsn.Pr
## Richness -0.623
## CmmntyMltpl 0.044 -0.611
## Fusn.Prtnrs -0.291 0.579 -0.673
## Rchnss:Fs.P 0.384 -0.723 0.582 -0.912
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(polypmod1, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Polyps
## Chisq Df Pr(>Chisq)
## Richness 0.2843 1 0.5939
## Community 0.3530 1 0.5524
## Fusion.Partners 0.5587 1 0.4548
## Richness:Fusion.Partners 0.0360 1 0.8494qqPlot(residuals(polypmod1))
## 389 393
## 73 74meangrowth <- ddply(final, c("Community"), summarise,
N = length(Polyps[!is.na(Polyps)]),
mean = mean(Polyps, na.rm=TRUE),
sd = sd(Polyps, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Polyps, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meangrowth## Community N mean sd se max lower upper
## 1 Individual 45 3.533333 1.778661 0.2651472 8 1.754672 5.311995
## 2 Multiple 112 3.357143 1.558800 0.1472928 7 1.798343 4.915943Plot the model outputs.
eff.polyps1 <- predictorEffect("Fusion.Partners", polypmod1, xlevels=4, rug=TRUE)
polyp.effects1<-plot(eff.polyps1,
lines=list(multiline=TRUE,
col=c("lightblue","mediumblue", "blue4", "darkblue"),
lty=1),
confint=list(style="bands", alpha=0.1),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("Total Polyp Expansion")),
legend.position="right",
xlab=expression(bold("Fusion Partners")),
main="",
lattice=list(key.args=list(space="right",
border=FALSE,
title=expression(bold("Genotypic \nRichness")),
cex=1,
cex.title=1)));polyp.effects1
Combine figures for Grow-Out Period.
grow_out_figs<-plot_grid(surv.effects1, surv.effects2b, labels = c("A", "B"), ncol=2, nrow=1, rel_heights= c(1,1), rel_widths = c(0.8,1), label_size = 20, label_y=1, align="h");grow_out_figs
ggsave(filename="Figures/grow_out_figs.png", plot=grow_out_figs, dpi=500, width=10, height=6, units="in")4. STRESS SURVIVORSHIP
Load data
Survivorship is measured as binary on each individual spat in the experiment.
#load data
stress<-read.csv("Data/Stress_Responses_PostGenotyping.csv", header=TRUE, sep=",", na.strings="NA")
stress$Diversity.Type<-relevel(stress$Diversity.Type, ref="Individual")
stress$Polyps<-as.numeric(stress$Polyps)Subset the final timepoint, “multiple” juvenile datasets for additional analyses.
final_stress<-stress[which(stress$Timepoint == "End"),]
multiple_stress<-stress[which(stress$Community == "Multiple"),]
multiple_final_stress<-multiple_stress[which(multiple_stress$Timepoint == "End"),]
d26_stress<-stress[which(stress$Timepoint == "S2D11"),]
d23_stress<-stress[which(stress$Timepoint == "S2D8"),]Analysis
Analyze with a logistic binomial regression model.
Check for effect of community type (“individual” or “multiple” juveniles) on survivorship and display means.
#with full time dataset
survmod3b<-glmer(Survivorship ~ Day * Treatment * Community + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial, link=logit)
summary(survmod3b)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula:
## Survivorship ~ Day * Treatment * Community + (1 | Parent.Site/Colony) +
## (1 | Treatment:Tank) + (1 | Slide)
## Data: stress
##
## AIC BIC logLik deviance df.resid
## 449.0 508.9 -212.5 425.0 1076
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -110.079 0.001 0.028 0.118 5.533
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 2.605e+00 1.614e+00
## Colony:Parent.Site (Intercept) 5.717e-01 7.561e-01
## Treatment:Tank (Intercept) 5.061e-09 7.114e-05
## Parent.Site (Intercept) 1.645e-08 1.282e-04
## Number of obs: 1088, groups:
## Slide, 90; Colony:Parent.Site, 20; Treatment:Tank, 9; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 6.7660321 1.7130404 3.950 7.82e-05 ***
## Day -0.0835370 0.0592967 -1.409 0.158896
## TreatmentHigh 0.8459104 1.9259267 0.439 0.660500
## CommunityMultiple 3.4775795 2.4771976 1.404 0.160368
## Day:TreatmentHigh -0.2886216 0.0823811 -3.503 0.000459 ***
## Day:CommunityMultiple -0.1855821 0.0878809 -2.112 0.034708 *
## TreatmentHigh:CommunityMultiple 2.8355215 3.2070137 0.884 0.376608
## Day:TreatmentHigh:CommunityMultiple 0.0001848 0.1226640 0.002 0.998798
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Day TrtmnH CmmntM Dy:TrH Dy:CmM TrH:CM
## Day -0.835
## TreatmntHgh -0.700 0.666
## CmmntyMltpl -0.586 0.538 0.484
## Dy:TrtmntHg 0.397 -0.633 -0.837 -0.348
## Dy:CmmntyMl 0.507 -0.652 -0.450 -0.917 0.471
## TrtmntHg:CM 0.412 -0.395 -0.600 -0.760 0.504 0.702
## Dy:TrtmH:CM -0.295 0.436 0.563 0.641 -0.648 -0.713 -0.924
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(survmod3b, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Survivorship
## Chisq Df Pr(>Chisq)
## Day 65.5353 1 5.708e-16 ***
## Treatment 37.2722 1 1.027e-09 ***
## Community 3.3205 1 0.068421 .
## Day:Treatment 21.1660 1 4.212e-06 ***
## Day:Community 9.0584 1 0.002615 **
## Treatment:Community 5.3796 1 0.020374 *
## Day:Treatment:Community 0.0000 1 0.998798
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1dispersion_glmer(survmod3b) #no overdispersion## [1] 0.5599041plot(residuals(survmod3b)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)meanstresssurv1 <- ddply(final_stress, c("Community", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv1## Community Treatment N mean sd se max lower
## 1 Individual Ambient 21 0.95238095 0.2182179 0.04761905 1 0.7341631
## 2 Individual High 24 0.00000000 0.0000000 0.00000000 0 0.0000000
## 3 Multiple Ambient 46 0.76086957 0.4312660 0.06358670 1 0.3296036
## 4 Multiple High 65 0.04615385 0.2114510 0.02622727 1 -0.1652972
## upper
## 1 1.1705988
## 2 0.0000000
## 3 1.1921355
## 4 0.2576049Analyze within “multiple” groups as survival was significantly different between single and multiple juvenile slides.
#with full time dataset
survmod3<-glmer(Survivorship ~ Day + Treatment + Fusion.Partners + Richness + Day:Treatment + Day:Fusion.Partners + Day:Richness + Day:Treatment:Fusion.Partners + Day:Treatment:Richness + Day:Fusion.Partners:Richness:Treatment + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, subset=c(Community=="Multiple"))
#fixed_form<-glm(Survivorship ~ Day * Treatment * Fusion.Partners *Richness, data=stress, family=binomial, subset=c(Community=="Multiple"))
#library("brglm2"); glm(fixed_form, data = stress, family = binomial, method="detect_separation")
summary(survmod3)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: binomial ( logit )
## Formula: Survivorship ~ Day + Treatment + Fusion.Partners + Richness +
## Day:Treatment + Day:Fusion.Partners + Day:Richness + Day:Treatment:Fusion.Partners +
## Day:Treatment:Richness + Day:Fusion.Partners:Richness:Treatment +
## (1 | Parent.Site/Colony) + (1 | Treatment:Tank) + (1 | Slide/Sample)
## Data: stress
## Subset: c(Community == "Multiple")
##
## AIC BIC logLik deviance df.resid
## 222.5 301.5 -94.2 188.5 756
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -0.6798 0.0000 0.0000 0.0000 0.5145
##
## Random effects:
## Groups Name Variance Std.Dev.
## Sample:Slide (Intercept) 940.14567 30.6618
## Slide (Intercept) 0.04049 0.2012
## Colony:Parent.Site (Intercept) 0.01752 0.1324
## Treatment:Tank (Intercept) 8.33415 2.8869
## Parent.Site (Intercept) 0.03312 0.1820
## Number of obs: 773, groups:
## Sample:Slide, 111; Slide, 45; Colony:Parent.Site, 20; Treatment:Tank, 8; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value
## (Intercept) -2.29443 39.97253 -0.057
## Day 0.55107 1.24309 0.443
## TreatmentHigh 111.46301 35.69621 3.123
## Fusion.Partners 100.44513 51.46355 1.952
## Richness 22.39796 34.13494 0.656
## Day:TreatmentHigh -4.98049 1.87839 -2.651
## Day:Fusion.Partners -2.73031 1.88657 -1.447
## Day:Richness -0.90135 1.21483 -0.742
## Day:TreatmentHigh:Fusion.Partners -1.42109 0.99729 -1.425
## Day:TreatmentHigh:Richness -0.02063 0.80677 -0.026
## Day:TreatmentAmbient:Fusion.Partners:Richness -0.09927 0.27596 -0.360
## Day:TreatmentHigh:Fusion.Partners:Richness 0.07110 0.22617 0.314
## Pr(>|z|)
## (Intercept) 0.95423
## Day 0.65754
## TreatmentHigh 0.00179 **
## Fusion.Partners 0.05097 .
## Richness 0.51172
## Day:TreatmentHigh 0.00801 **
## Day:Fusion.Partners 0.14783
## Day:Richness 0.45811
## Day:TreatmentHigh:Fusion.Partners 0.15417
## Day:TreatmentHigh:Richness 0.97959
## Day:TreatmentAmbient:Fusion.Partners:Richness 0.71905
## Day:TreatmentHigh:Fusion.Partners:Richness 0.75326
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Day TrtmnH Fsn.Pr Rchnss Dy:TrH Dy:F.P Dy:Rch Dy:TH:F.P
## Day -0.633
## TreatmntHgh 0.030 -0.317
## Fusn.Prtnrs 0.516 -0.509 0.005
## Richness -0.917 0.826 -0.219 -0.590
## Dy:TrtmntHg -0.392 0.059 -0.623 -0.043 0.331
## Dy:Fsn.Prtn -0.611 0.393 0.087 -0.936 0.590 0.133
## Day:Richnss 0.627 -0.990 0.337 0.519 -0.838 -0.074 -0.404
## Dy:TrtH:F.P 0.061 0.336 -0.155 -0.324 0.137 -0.338 0.071 -0.337
## Dy:TrtmnH:R 0.478 0.213 -0.071 0.098 -0.287 -0.680 -0.275 -0.210 0.453
## Dy:TA:F.P:R 0.424 0.315 -0.245 0.088 -0.181 -0.371 -0.391 -0.313 0.552
## Dy:TH:F.P:R 0.162 -0.148 -0.066 0.147 -0.178 0.307 -0.142 0.151 -0.652
## D:TH:R D:TA:F
## Day
## TreatmntHgh
## Fusn.Prtnrs
## Richness
## Dy:TrtmntHg
## Dy:Fsn.Prtn
## Day:Richnss
## Dy:TrtH:F.P
## Dy:TrtmnH:R
## Dy:TA:F.P:R 0.752
## Dy:TH:F.P:R -0.255 0.031
## convergence code: 0
## unable to evaluate scaled gradient
## Model failed to converge: degenerate Hessian with 6 negative eigenvaluesAnova(survmod3, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Survivorship
## Chisq Df Pr(>Chisq)
## Day 8.3067 1 0.00395 **
## Treatment 5.5318 1 0.01867 *
## Fusion.Partners 0.6676 1 0.41390
## Richness 0.2974 1 0.58553
## Day:Treatment 16.8143 1 4.122e-05 ***
## Day:Fusion.Partners 0.8008 1 0.37084
## Day:Richness 1.0944 1 0.29549
## Day:Treatment:Fusion.Partners 4.1122 1 0.04257 *
## Day:Treatment:Richness 0.3144 1 0.57499
## Day:Treatment:Fusion.Partners:Richness 0.2355 2 0.88894
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1dispersion_glmer(survmod3) #no overdispersion## [1] 0.2731621plot(residuals(survmod3)) + abline(h=0, lty=2) #Dispersed randomly
## integer(0)meanstresssurv2 <- ddply(final_stress, c("Fusion.Partners", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv2## Fusion.Partners Treatment N mean sd se max lower
## 1 0 Ambient 30 0.83333333 0.3790490 0.06920457 1 0.4542843
## 2 0 High 50 0.04000000 0.1979487 0.02799417 1 -0.1579487
## 3 1 Ambient 16 0.81250000 0.4031129 0.10077822 1 0.4093871
## 4 1 High 15 0.06666667 0.2581989 0.06666667 1 -0.1915322
## 5 2 Ambient 9 0.77777778 0.4409586 0.14698618 1 0.3368192
## 6 2 High 16 0.00000000 0.0000000 0.00000000 0 0.0000000
## 7 3 Ambient 12 0.83333333 0.3892495 0.11236664 1 0.4440839
## 8 3 High 8 0.00000000 0.0000000 0.00000000 0 0.0000000
## upper
## 1 1.2123824
## 2 0.2379487
## 3 1.2156129
## 4 0.3248656
## 5 1.2187363
## 6 0.0000000
## 7 1.2225828
## 8 0.0000000Plot effect of number of juveniles over time.
stress$group<-paste(stress$Treatment, "-", stress$Community)
survmod3b2<-glmer(Survivorship ~ group * Day + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial, link=logit)
eff.surv3b2 <- predictorEffect(c("Day"), survmod3b2)
surv.effects3b2<-plot(eff.surv3b2,
lines=list(multiline=TRUE,
col=c("blue", "blue", "red", "red"),
lty=c(1, 2, 1, 2)),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1.1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="top",
xlab=expression(bold("Day")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Temperature - Number of Juveniles")),
cex=1,
cex.title=1)));surv.effects3b2
Examine curve estimates to identify 50% survivorship.
summary(effect(c("group", "Day"), xlevels=50, survmod3b2))##
## group*Day effect
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9989505 0.9988891 0.9988237 0.9987553 0.9986828
## Ambient - Multiple 0.9999639 0.9999570 0.9999487 0.9999389 0.9999273
## High - Individual 0.9995124 0.9993782 0.9992059 0.9989885 0.9987118
## High - Multiple 0.9999991 0.9999986 0.9999980 0.9999972 0.9999960
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9986050 0.9985238 0.9984379 0.9983470 0.9982494
## Ambient - Multiple 0.9999132 0.9998966 0.9998769 0.9998534 0.9998250
## High - Individual 0.9983533 0.9979030 0.9973299 0.9966008 0.9956574
## High - Multiple 0.9999942 0.9999916 0.9999880 0.9999828 0.9999752
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9981475 0.9980398 0.9979240 0.9978033 0.9976756
## Ambient - Multiple 0.9997916 0.9997518 0.9997037 0.9996472 0.9995799
## High - Individual 0.9944741 0.9929706 0.9910287 0.9885985 0.9855198
## High - Multiple 0.9999644 0.9999489 0.9999263 0.9998942 0.9998482
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9975384 0.9974066 0.9972439 0.9970711 0.9969144
## Ambient - Multiple 0.9994984 0.9994107 0.9992889 0.9991419 0.9989919
## High - Individual 0.9815577 0.9770442 0.9704068 0.9619250 0.9528374
## High - Multiple 0.9997810 0.9996944 0.9995493 0.9993353 0.9990726
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9967210 0.9965456 0.9963292 0.9961329 0.9958907
## Ambient - Multiple 0.9987836 0.9985711 0.9982760 0.9979751 0.9975572
## High - Individual 0.9396333 0.9256365 0.9055717 0.8846468 0.8552516
## High - Multiple 0.9986325 0.9980927 0.9971889 0.9960815 0.9942302
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9956711 0.9954001 0.9951544 0.9948512 0.9945764
## Ambient - Multiple 0.9971311 0.9965396 0.9959368 0.9951002 0.9942481
## High - Individual 0.8253269 0.7844978 0.7443186 0.6916290 0.6420357
## High - Multiple 0.9919665 0.9881940 0.9836013 0.9759955 0.9668167
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9942373 0.9939299 0.9935506 0.9931477 0.9927826
## Ambient - Multiple 0.9930663 0.9918634 0.9901965 0.9881921 0.9861554
## High - Individual 0.5801574 0.5249503 0.4598598 0.3961105 0.3440645
## High - Multiple 0.9518077 0.9340054 0.9056046 0.8667256 0.8233276
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9923322 0.9919240 0.9914204 0.9909641 0.9904013
## Ambient - Multiple 0.9833386 0.9804813 0.9765376 0.9725465 0.9670535
## High - Individual 0.2878140 0.2442424 0.1993517 0.1660498 0.1330014
## High - Multiple 0.7595582 0.6935999 0.6054450 0.5237196 0.4270620
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.9898913 0.98926241 0.98869264 0.98799004 0.9873536
## Ambient - Multiple 0.9615127 0.95391658 0.94628912 0.93588877 0.9255105
## High - Individual 0.1092707 0.08635276 0.07027067 0.05502712 0.0444949
## High - Multiple 0.3481675 0.26582641 0.20600815 0.14957280 0.1119266
## Day
## group 29.4 30 30.7 31.3 32
## Ambient - Individual 0.98656895 0.98585830 0.98498226 0.98418900 0.98321132
## Ambient - Multiple 0.91146307 0.89756352 0.87893523 0.86070948 0.83659930
## High - Individual 0.03463450 0.02789024 0.02162627 0.01736948 0.01343575
## High - Multiple 0.07870967 0.05768922 0.03984634 0.02887951 0.01976041
##
## Lower 95 Percent Confidence Limits
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9688623 0.9690922 0.9693058 0.9694989 0.9696732
## Ambient - Multiple 0.9981187 0.9979334 0.9977283 0.9975049 0.9972591
## High - Individual 0.9920845 0.9906655 0.9889803 0.9870133 0.9846947
## High - Multiple 0.9999701 0.9999606 0.9999479 0.9999313 0.9999094
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9698303 0.9699644 0.9700772 0.9701676 0.9702353
## Ambient - Multiple 0.9969842 0.9966859 0.9963574 0.9959955 0.9955902
## High - Individual 0.9819168 0.9786896 0.9748889 0.9704148 0.9650636
## High - Multiple 0.9998801 0.9998418 0.9997913 0.9997246 0.9996350
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9702769 0.9702921 0.9702789 0.9702358 0.9701607
## Ambient - Multiple 0.9951498 0.9946641 0.9941193 0.9935266 0.9928718
## High - Individual 0.9588613 0.9515783 0.9428952 0.9328710 0.9211579
## High - Multiple 0.9995181 0.9993637 0.9991558 0.9988845 0.9985255
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9700495 0.9699156 0.9697168 0.9694682 0.9692122
## Ambient - Multiple 0.9921364 0.9913992 0.9904473 0.9893846 0.9883749
## High - Individual 0.9072779 0.8926825 0.8729571 0.8499648 0.8273689
## High - Multiple 0.9980419 0.9974650 0.9965725 0.9953636 0.9939910
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9688590 0.9685055 0.9680282 0.9675585 0.9669328
## Ambient - Multiple 0.9870676 0.9858233 0.9842088 0.9826685 0.9806649
## High - Individual 0.7973226 0.7681915 0.7300762 0.6938098 0.6473968
## High - Multiple 0.9918652 0.9894505 0.9857101 0.9814630 0.9748907
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9663236 0.9655190 0.9647409 0.9637193 0.9627359
## Ambient - Multiple 0.9787481 0.9762472 0.9738466 0.9707029 0.9676733
## High - Individual 0.6043419 0.5508357 0.5028190 0.4453401 0.3958471
## High - Multiple 0.9674434 0.9559596 0.9430171 0.9232269 0.9011890
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9614496 0.9602154 0.9586051 0.9567911 0.9550561
## Ambient - Multiple 0.9636884 0.9598297 0.9547280 0.9488733 0.9431543
## High - Individual 0.3392213 0.2927669 0.2422482 0.1967494 0.1622952
## High - Multiple 0.8680761 0.8320749 0.7797593 0.7153722 0.6506739
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9527987 0.9506428 0.94784145 0.94516956 0.9417021
## Ambient - Multiple 0.9355228 0.9280213 0.91794801 0.90798712 0.8945388
## High - Individual 0.1276564 0.1026581 0.07860727 0.06193355 0.0464441
## High - Multiple 0.5657819 0.4877617 0.39554313 0.31998675 0.2407995
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.93839918 0.93411856 0.93004678 0.92477792 0.91977455
## Ambient - Multiple 0.88118186 0.86309627 0.84511815 0.82082232 0.79679637
## High - Individual 0.03603222 0.02660821 0.02041373 0.01490956 0.01134799
## High - Multiple 0.18322882 0.12923146 0.09371102 0.06306962 0.04427288
## Day
## group 29.4 30 30.7 31.3
## Ambient - Individual 0.913312333 0.907188427 0.899297391 0.891838810
## Ambient - Multiple 0.764632262 0.733264999 0.692078154 0.652889335
## High - Individual 0.008223676 0.006223925 0.004485222 0.003380753
## High - Multiple 0.028915333 0.019893620 0.012757874 0.008671295
## Day
## group 32
## Ambient - Individual 0.882256056
## Ambient - Multiple 0.602984743
## High - Individual 0.002426434
## High - Multiple 0.005499545
##
## Upper 95 Percent Confidence Limits
## Day
## group 0 0.653 1.31 1.96 2.61
## Ambient - Individual 0.9999657 0.9999612 0.9999562 0.9999506 0.9999444
## Ambient - Multiple 0.9999993 0.9999991 0.9999988 0.9999985 0.9999981
## High - Individual 0.9999702 0.9999589 0.9999433 0.9999221 0.9998930
## High - Multiple 1.0000000 1.0000000 0.9999999 0.9999999 0.9999998
## Day
## group 3.27 3.92 4.57 5.22 5.88
## Ambient - Individual 0.9999373 0.9999294 0.9999207 0.9999109 0.9998998
## Ambient - Multiple 0.9999975 0.9999968 0.9999959 0.9999947 0.9999931
## High - Individual 0.9998523 0.9997972 0.9997218 0.9996186 0.9994748
## High - Multiple 0.9999997 0.9999996 0.9999993 0.9999989 0.9999983
## Day
## group 6.53 7.18 7.84 8.49 9.14
## Ambient - Individual 0.9998876 0.9998740 0.9998587 0.9998420 0.9998236
## Ambient - Multiple 0.9999911 0.9999885 0.9999851 0.9999809 0.9999754
## High - Individual 0.9992809 0.9990161 0.9986487 0.9981550 0.9974840
## High - Multiple 0.9999974 0.9999959 0.9999936 0.9999900 0.9999844
## Day
## group 9.8 10.4 11.1 11.8 12.4
## Ambient - Individual 0.9998028 0.9997821 0.9997555 0.9997261 0.9996985
## Ambient - Multiple 0.9999682 0.9999599 0.9999475 0.9999313 0.9999134
## High - Individual 0.9965576 0.9954292 0.9936503 0.9912023 0.9883944
## High - Multiple 0.9999755 0.9999632 0.9999409 0.9999050 0.9998575
## Day
## group 13.1 13.7 14.4 15 15.7
## Ambient - Individual 0.9996634 0.9996306 0.9995892 0.9995507 0.9995024
## Ambient - Multiple 0.9998868 0.9998576 0.9998142 0.9997666 0.9996959
## High - Individual 0.9840224 0.9790595 0.9714308 0.9629020 0.9500353
## High - Multiple 0.9997714 0.9996576 0.9994521 0.9991813 0.9986941
## Day
## group 16.3 17 17.6 18.3 18.9
## Ambient - Individual 0.9994579 0.9994024 0.9993517 0.9992890 0.9992323
## Ambient - Multiple 0.9996189 0.9995047 0.9993806 0.9991973 0.9989992
## High - Individual 0.9359643 0.9152970 0.8933854 0.8623573 0.8307878
## High - Multiple 0.9980548 0.9969114 0.9954211 0.9927783 0.9893703
## Day
## group 19.6 20.2 20.9 21.6 22.2
## Ambient - Individual 0.9991628 0.9991006 0.9990252 0.9989470 0.9988782
## Ambient - Multiple 0.9987079 0.9983946 0.9979371 0.9973571 0.9967406
## High - Individual 0.7881158 0.7468255 0.6939347 0.6372250 0.5868072
## High - Multiple 0.9834108 0.9758592 0.9629580 0.9439063 0.9210072
## Day
## group 22.9 23.5 24.2 24.8 25.5
## Ambient - Individual 0.9987962 0.9987249 0.9986410 0.9985688 0.9984850
## Ambient - Multiple 0.9958518 0.9949165 0.9935835 0.9921980 0.9902508
## High - Individual 0.5274205 0.4772442 0.4208513 0.3751904 0.3257642
## High - Multiple 0.8845106 0.8432974 0.7825316 0.7198529 0.6365925
## Day
## group 26.1 26.8 27.4 28.1 28.7
## Ambient - Individual 0.9984139 0.9983324 0.9982640 0.9981867 0.9981227
## Ambient - Multiple 0.9882570 0.9855000 0.9827247 0.9789551 0.9752286
## High - Individual 0.2870447 0.2463006 0.2151503 0.1830334 0.1589008
## High - Multiple 0.5598160 0.4690315 0.3943239 0.3148500 0.2553418
## Day
## group 29.4 30 30.7 31.3 32
## Ambient - Individual 0.9980511 0.9979928 0.99792836 0.99787649 0.99782004
## Ambient - Multiple 0.9702584 0.9654314 0.95910208 0.95305189 0.94523397
## High - Individual 0.1343732 0.1161636 0.09783686 0.08434176 0.07084932
## High - Multiple 0.1968694 0.1558722 0.11759973 0.09182032 0.06845565Create plot for effect of temperature on survival in multiple juvenile groups across number of fusion partners.
stress$group<-paste(stress$Treatment, "-", stress$Fusion.Partners)
survmod3a<-glmer(Survivorship ~ Day * group + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=stress, family=binomial)
summary(effect(c("Day", "group"), xlevels=50, survmod3a))##
## Day*group effect
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0 High - 1
## 0 0.9996515 0.9999990 0.9999988 0.9999991 0.99974119 0.9999999995
## 0.653 0.9996149 0.9999987 0.9999985 0.9999989 0.99967271 0.9999999990
## 1.31 0.9995742 0.9999983 0.9999981 0.9999986 0.99958552 0.9999999982
## 1.96 0.9995297 0.9999979 0.9999975 0.9999983 0.99947644 0.9999999967
## 2.61 0.9994806 0.9999973 0.9999969 0.9999978 0.99933866 0.9999999939
## 3.27 0.9994255 0.9999965 0.9999960 0.9999973 0.99916166 0.9999999888
## 3.92 0.9993654 0.9999956 0.9999949 0.9999966 0.99894113 0.9999999796
## 4.57 0.9992992 0.9999944 0.9999935 0.9999957 0.99866268 0.9999999629
## 5.22 0.9992260 0.9999929 0.9999917 0.9999945 0.99831112 0.9999999324
## 5.88 0.9991438 0.9999909 0.9999895 0.9999931 0.99785968 0.9999998757
## 6.53 0.9990544 0.9999884 0.9999866 0.9999912 0.99729760 0.9999997737
## 7.18 0.9989557 0.9999853 0.9999830 0.9999890 0.99658842 0.9999995879
## 7.84 0.9988449 0.9999812 0.9999782 0.9999860 0.99567848 0.9999992427
## 8.49 0.9987243 0.9999762 0.9999723 0.9999824 0.99454672 0.9999986208
## 9.14 0.9985911 0.9999697 0.9999648 0.9999778 0.99312061 0.9999974885
## 9.8 0.9984417 0.9999614 0.9999551 0.9999718 0.99129383 0.9999953841
## 10.4 0.9982922 0.9999518 0.9999439 0.9999651 0.98921987 0.9999919730
## 11.1 0.9980996 0.9999377 0.9999273 0.9999552 0.98617675 0.9999846928
## 11.8 0.9978853 0.9999193 0.9999058 0.9999425 0.98228997 0.9999708098
## 12.4 0.9976825 0.9998994 0.9998824 0.9999288 0.97811848 0.9999492398
## 13.1 0.9974213 0.9998697 0.9998477 0.9999085 0.97203016 0.9999032059
## 13.7 0.9971741 0.9998375 0.9998098 0.9998867 0.96552684 0.9998316884
## 14.4 0.9968558 0.9997897 0.9997536 0.9998545 0.95609229 0.9996790825
## 15 0.9965547 0.9997377 0.9996924 0.9998197 0.94608877 0.9994420611
## 15.7 0.9961668 0.9996605 0.9996015 0.9997686 0.93171070 0.9989365612
## 16.3 0.9957999 0.9995766 0.9995024 0.9997133 0.91663663 0.9981521464
## 17 0.9953275 0.9994521 0.9993554 0.9996319 0.89527321 0.9964820804
## 17.6 0.9948807 0.9993166 0.9991953 0.9995440 0.87325018 0.9938982477
## 18.3 0.9943055 0.9991157 0.9989577 0.9994146 0.84267645 0.9884281352
## 18.9 0.9937616 0.9988971 0.9986989 0.9992748 0.81191650 0.9800473070
## 19.6 0.9930615 0.9985730 0.9983149 0.9990690 0.77043751 0.9626269868
## 20.2 0.9923997 0.9982206 0.9978969 0.9988469 0.73007757 0.9367555273
## 20.9 0.9915481 0.9976982 0.9972768 0.9985199 0.67771441 0.8859370088
## 21.6 0.9906020 0.9970230 0.9964745 0.9981004 0.62047378 0.8028772805
## 22.2 0.9897080 0.9962891 0.9956020 0.9976476 0.56851442 0.7007921559
## 22.9 0.9885583 0.9952024 0.9943092 0.9969815 0.50601579 0.5512081429
## 23.5 0.9874724 0.9940225 0.9929045 0.9962630 0.45222054 0.4139283323
## 24.2 0.9860764 0.9922772 0.9908260 0.9952069 0.39092373 0.2702669052
## 24.8 0.9847586 0.9903848 0.9885716 0.9940685 0.34091977 0.1755820234
## 25.5 0.9830654 0.9875907 0.9852427 0.9923969 0.28681031 0.1004629759
## 26.1 0.9814679 0.9845677 0.9816416 0.9905973 0.24477183 0.0603473376
## 26.8 0.9794168 0.9801171 0.9763424 0.9879594 0.20126260 0.0325806638
## 27.4 0.9774830 0.9753189 0.9706342 0.9851252 0.16879538 0.0189984697
## 28.1 0.9750021 0.9682869 0.9622797 0.9809815 0.13635299 0.0100534491
## 28.7 0.9726651 0.9607483 0.9533406 0.9765440 0.11287708 0.0058060161
## 29.4 0.9696698 0.9497787 0.9403684 0.9700824 0.09001823 0.0030530475
## 30 0.9668513 0.9381216 0.9266326 0.9631971 0.07383774 0.0017579282
## 30.7 0.9632430 0.9213463 0.9069597 0.9532348 0.05836465 0.0009226114
## 31.3 0.9598521 0.9037583 0.8864575 0.9427013 0.04757632 0.0005307535
## 32 0.9555170 0.8788699 0.8576668 0.9276079 0.03738431 0.0002783921
## group
## Day High - 2 High - 3
## 0 1.000000e+00 1.000000e+00
## 0.653 1.000000e+00 1.000000e+00
## 1.31 1.000000e+00 1.000000e+00
## 1.96 1.000000e+00 1.000000e+00
## 2.61 1.000000e+00 1.000000e+00
## 3.27 1.000000e+00 1.000000e+00
## 3.92 1.000000e+00 1.000000e+00
## 4.57 1.000000e+00 1.000000e+00
## 5.22 1.000000e+00 1.000000e+00
## 5.88 1.000000e+00 1.000000e+00
## 6.53 1.000000e+00 1.000000e+00
## 7.18 1.000000e+00 1.000000e+00
## 7.84 1.000000e+00 1.000000e+00
## 8.49 1.000000e+00 1.000000e+00
## 9.14 1.000000e+00 1.000000e+00
## 9.8 1.000000e+00 1.000000e+00
## 10.4 1.000000e+00 1.000000e+00
## 11.1 1.000000e+00 1.000000e+00
## 11.8 1.000000e+00 1.000000e+00
## 12.4 1.000000e+00 1.000000e+00
## 13.1 1.000000e+00 1.000000e+00
## 13.7 1.000000e+00 1.000000e+00
## 14.4 1.000000e+00 1.000000e+00
## 15 1.000000e+00 9.999999e-01
## 15.7 1.000000e+00 9.999997e-01
## 16.3 9.999999e-01 9.999993e-01
## 17 9.999996e-01 9.999978e-01
## 17.6 9.999988e-01 9.999941e-01
## 18.3 9.999956e-01 9.999816e-01
## 18.9 9.999862e-01 9.999515e-01
## 19.6 9.999476e-01 9.998494e-01
## 20.2 9.998352e-01 9.996025e-01
## 20.9 9.993730e-01 9.987669e-01
## 21.6 9.976180e-01 9.961820e-01
## 22.2 9.925485e-01 9.899812e-01
## 22.9 9.722222e-01 9.695402e-01
## 23.5 9.175687e-01 9.233977e-01
## 24.2 7.452137e-01 7.952102e-01
## 24.8 4.819230e-01 5.952330e-01
## 25.5 1.964146e-01 3.214383e-01
## 26.1 7.212871e-02 1.521096e-01
## 26.8 2.001693e-02 5.463164e-02
## 27.4 6.454244e-03 2.141655e-02
## 28.1 1.704022e-03 7.000467e-03
## 28.7 5.425739e-04 2.662738e-03
## 29.4 1.426235e-04 8.592894e-04
## 30 4.536407e-05 3.255965e-04
## 30.7 1.192024e-05 1.049064e-04
## 31.3 3.791117e-06 3.973184e-05
## 32 9.961548e-07 1.279902e-05
##
## Lower 95 Percent Confidence Limits
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0 High - 1
## 0 0.9912555 0.9502673 0.7165756 0.8366981 0.99715478 9.996945e-01
## 0.653 0.9908775 0.9496746 0.7205324 0.8392318 0.99660966 9.996073e-01
## 1.31 0.9904778 0.9490632 0.7244073 0.8417018 0.99595527 9.994944e-01
## 1.96 0.9900620 0.9484425 0.7281329 0.8440663 0.99518309 9.993506e-01
## 2.61 0.9896245 0.9478050 0.7317479 0.8463507 0.99426285 9.991659e-01
## 3.27 0.9891567 0.9471393 0.7353016 0.8485870 0.99314769 9.989242e-01
## 3.92 0.9886710 0.9464645 0.7386825 0.8507058 0.99183702 9.986176e-01
## 4.57 0.9881589 0.9457692 0.7419412 0.8527396 0.99027491 9.982233e-01
## 5.22 0.9876185 0.9450521 0.7450728 0.8546865 0.98841323 9.977162e-01
## 5.88 0.9870386 0.9442998 0.7481175 0.8565718 0.98615752 9.970525e-01
## 6.53 0.9864345 0.9435333 0.7509773 0.8583356 0.98350730 9.962100e-01
## 7.18 0.9857952 0.9427393 0.7536930 0.8600040 0.98035072 9.951259e-01
## 7.84 0.9851071 0.9419028 0.7562959 0.8615967 0.97652946 9.937067e-01
## 8.49 0.9843881 0.9410463 0.7586989 0.8630612 0.97204557 9.919051e-01
## 9.14 0.9836246 0.9401546 0.7609340 0.8644176 0.96671406 9.895876e-01
## 9.8 0.9828000 0.9392098 0.7630212 0.8656786 0.96027466 9.865554e-01
## 10.4 0.9820036 0.9383131 0.7647492 0.8667176 0.95336712 9.830417e-01
## 11.1 0.9810131 0.9372170 0.7665473 0.8677926 0.94381724 9.777736e-01
## 11.8 0.9799506 0.9360613 0.7680929 0.8687096 0.93238126 9.708875e-01
## 12.4 0.9789770 0.9350180 0.7692005 0.8693606 0.92082640 9.633385e-01
## 13.1 0.9777606 0.9337320 0.7702175 0.8699494 0.90497463 9.520862e-01
## 13.7 0.9766426 0.9325641 0.7708325 0.8702956 0.88906659 9.398238e-01
## 14.4 0.9752412 0.9311152 0.7712210 0.8704969 0.86742606 9.216944e-01
## 15 0.9739487 0.9297899 0.7712440 0.8704788 0.84592741 9.021463e-01
## 15.7 0.9723229 0.9281332 0.7708687 0.8702106 0.81704140 8.736528e-01
## 16.3 0.9708179 0.9266054 0.7701635 0.8697450 0.78876238 8.434788e-01
## 17 0.9689176 0.9246780 0.7688371 0.8688921 0.75143047 8.005154e-01
## 17.6 0.9671517 0.9228830 0.7672141 0.8678612 0.71562740 7.563176e-01
## 18.3 0.9649131 0.9205943 0.7646734 0.8662585 0.66949610 6.956332e-01
## 18.9 0.9628244 0.9184380 0.7618631 0.8644918 0.62646612 6.358615e-01
## 19.6 0.9601659 0.9156535 0.7577311 0.8618979 0.57276614 5.579737e-01
## 20.2 0.9576754 0.9129939 0.7533448 0.8591434 0.52443219 4.857392e-01
## 20.9 0.9544923 0.9095072 0.7470733 0.8551976 0.46646243 3.979160e-01
## 21.6 0.9509737 0.9055171 0.7393365 0.8503119 0.40820365 3.102373e-01
## 22.2 0.9476572 0.9016048 0.7313237 0.8452262 0.35924119 2.388110e-01
## 22.9 0.9433935 0.8963268 0.7200574 0.8380248 0.30462141 1.642022e-01
## 23.5 0.9393606 0.8910428 0.7084547 0.8305407 0.26092759 1.109629e-01
## 24.2 0.9341585 0.8837522 0.6922050 0.8199334 0.21452202 6.395642e-02
## 24.8 0.9292231 0.8762793 0.6755256 0.8088816 0.17914515 3.672233e-02
## 25.5 0.9228400 0.8657091 0.6522578 0.7931644 0.14324933 1.760807e-02
## 26.1 0.9167706 0.8546003 0.6285129 0.7767364 0.11703505 8.815467e-03
## 26.8 0.9089073 0.8384952 0.5957030 0.7533319 0.09145063 3.725555e-03
## 27.4 0.9014216 0.8211872 0.5627122 0.7288949 0.07341441 1.721273e-03
## 28.1 0.8917183 0.7956362 0.5181564 0.6942884 0.05634446 6.798659e-04
## 28.7 0.8824824 0.7678786 0.4747835 0.6586120 0.04463258 3.013223e-04
## 29.4 0.8705217 0.7269205 0.4187943 0.6092226 0.03379979 1.148955e-04
## 30 0.8591572 0.6831003 0.3673832 0.5600263 0.02651346 4.980172e-05
## 30.7 0.8444797 0.6207078 0.3057136 0.4952240 0.01988446 1.862178e-05
## 31.3 0.8305846 0.5577061 0.2537627 0.4347508 0.01548829 7.968705e-06
## 32 0.8127241 0.4752420 0.1972006 0.3613996 0.01153507 2.945008e-06
## group
## Day High - 2 High - 3
## 0 1.000000e+00 9.997782e-01
## 0.653 1.000000e+00 9.997244e-01
## 1.31 1.000000e+00 9.996571e-01
## 1.96 1.000000e+00 9.995742e-01
## 2.61 1.000000e+00 9.994711e-01
## 3.27 1.000000e+00 9.993408e-01
## 3.92 9.999999e-01 9.991808e-01
## 4.57 9.999999e-01 9.989817e-01
## 5.22 9.999998e-01 9.987339e-01
## 5.88 9.999996e-01 9.984199e-01
## 6.53 9.999994e-01 9.980340e-01
## 7.18 9.999990e-01 9.975529e-01
## 7.84 9.999983e-01 9.969425e-01
## 8.49 9.999971e-01 9.961911e-01
## 9.14 9.999951e-01 9.952527e-01
## 9.8 9.999916e-01 9.940600e-01
## 10.4 9.999864e-01 9.927139e-01
## 11.1 9.999759e-01 9.907475e-01
## 11.8 9.999575e-01 9.882417e-01
## 12.4 9.999308e-01 9.855511e-01
## 13.1 9.998775e-01 9.816101e-01
## 13.7 9.998002e-01 9.773708e-01
## 14.4 9.996459e-01 9.711491e-01
## 15 9.994212e-01 9.644429e-01
## 15.7 9.989713e-01 9.545798e-01
## 16.3 9.983137e-01 9.439250e-01
## 17 9.969922e-01 9.282183e-01
## 17.6 9.950516e-01 9.112107e-01
## 18.3 9.911338e-01 8.860763e-01
## 18.9 9.853546e-01 8.587915e-01
## 19.6 9.736495e-01 8.183675e-01
## 20.2 9.563854e-01 7.743915e-01
## 20.9 9.216908e-01 7.091908e-01
## 21.6 8.607833e-01 6.251733e-01
## 22.2 7.761912e-01 5.351789e-01
## 22.9 6.268973e-01 4.081918e-01
## 23.5 4.568730e-01 2.863192e-01
## 24.2 2.459640e-01 1.516768e-01
## 24.8 1.076049e-01 6.851122e-02
## 25.5 2.816542e-02 1.995146e-02
## 26.1 6.914230e-03 5.580647e-03
## 26.8 1.119779e-03 1.061058e-03
## 27.4 2.150945e-04 2.319496e-04
## 28.1 2.957691e-05 3.659250e-05
## 28.7 5.233405e-06 7.208943e-06
## 29.4 6.779846e-07 1.048428e-06
## 30 1.160429e-07 1.968729e-07
## 30.7 1.464006e-08 2.751402e-08
## 31.3 2.465852e-09 5.039492e-09
## 32 3.068533e-10 6.891137e-10
##
## Upper 95 Percent Confidence Limits
## group
## Day Ambient - 0 Ambient - 1 Ambient - 2 Ambient - 3 High - 0 High - 1
## 0 0.9999862 1.0000000 1.0000000 1.0000000 0.9999765 1.00000000
## 0.653 0.9999839 1.0000000 1.0000000 1.0000000 0.9999685 1.00000000
## 1.31 0.9999811 1.0000000 1.0000000 1.0000000 0.9999577 1.00000000
## 1.96 0.9999779 1.0000000 1.0000000 1.0000000 0.9999433 1.00000000
## 2.61 0.9999742 1.0000000 1.0000000 1.0000000 0.9999241 1.00000000
## 3.27 0.9999699 1.0000000 1.0000000 1.0000000 0.9998980 1.00000000
## 3.92 0.9999648 1.0000000 1.0000000 1.0000000 0.9998635 1.00000000
## 4.57 0.9999590 1.0000000 1.0000000 1.0000000 0.9998174 1.00000000
## 5.22 0.9999521 1.0000000 1.0000000 1.0000000 0.9997559 1.00000000
## 5.88 0.9999441 1.0000000 1.0000000 1.0000000 0.9996724 1.00000000
## 6.53 0.9999349 1.0000000 1.0000000 1.0000000 0.9995623 1.00000000
## 7.18 0.9999242 1.0000000 1.0000000 1.0000000 0.9994157 1.00000000
## 7.84 0.9999115 1.0000000 1.0000000 1.0000000 0.9992168 1.00000000
## 8.49 0.9998971 1.0000000 1.0000000 1.0000000 0.9989556 1.00000000
## 9.14 0.9998805 1.0000000 1.0000000 1.0000000 0.9986084 1.00000000
## 9.8 0.9998608 1.0000000 1.0000000 1.0000000 0.9981389 1.00000000
## 10.4 0.9998403 1.0000000 1.0000000 1.0000000 0.9975780 1.00000000
## 11.1 0.9998127 0.9999999 1.0000000 1.0000000 0.9967102 0.99999999
## 11.8 0.9997805 0.9999999 1.0000000 1.0000000 0.9955379 0.99999997
## 12.4 0.9997488 0.9999999 1.0000000 1.0000000 0.9942131 0.99999993
## 13.1 0.9997062 0.9999998 0.9999999 0.9999999 0.9921764 0.99999981
## 13.7 0.9996643 0.9999996 0.9999999 0.9999999 0.9898867 0.99999956
## 14.4 0.9996083 0.9999994 0.9999998 0.9999999 0.9863886 0.99999879
## 15 0.9995533 0.9999991 0.9999997 0.9999998 0.9824843 0.99999713
## 15.7 0.9994801 0.9999985 0.9999995 0.9999996 0.9765718 0.99999216
## 16.3 0.9994085 0.9999977 0.9999992 0.9999995 0.9700414 0.99998153
## 17 0.9993135 0.9999963 0.9999986 0.9999991 0.9602770 0.99994999
## 17.6 0.9992210 0.9999944 0.9999979 0.9999986 0.9496521 0.99988304
## 18.3 0.9990988 0.9999909 0.9999965 0.9999978 0.9340511 0.99968684
## 18.9 0.9989804 0.9999863 0.9999946 0.9999966 0.9174305 0.99927675
## 19.6 0.9988247 0.9999778 0.9999911 0.9999946 0.8936350 0.99810093
## 20.2 0.9986746 0.9999667 0.9999864 0.9999919 0.8690091 0.99571307
## 20.9 0.9984784 0.9999465 0.9999780 0.9999870 0.8349231 0.98916357
## 21.6 0.9982572 0.9999146 0.9999645 0.9999794 0.7948669 0.97360319
## 22.2 0.9980460 0.9998729 0.9999469 0.9999696 0.7558840 0.94590269
## 22.9 0.9977724 0.9997991 0.9999157 0.9999526 0.7054766 0.88476971
## 23.5 0.9975130 0.9997044 0.9998759 0.9999310 0.6587542 0.79986403
## 24.2 0.9971792 0.9995397 0.9998072 0.9998944 0.6013319 0.66750510
## 24.8 0.9968649 0.9993329 0.9997218 0.9998493 0.5507629 0.54334440
## 25.5 0.9964634 0.9989832 0.9995794 0.9997750 0.4916777 0.41034263
## 26.1 0.9960882 0.9985581 0.9994086 0.9996867 0.4421170 0.31682661
## 26.8 0.9956125 0.9978680 0.9991356 0.9995466 0.3867976 0.23271938
## 27.4 0.9951711 0.9970677 0.9988235 0.9993874 0.3423158 0.17865853
## 28.1 0.9946157 0.9958412 0.9983504 0.9991471 0.2945147 0.13164000
## 28.7 0.9941042 0.9945086 0.9978393 0.9988882 0.2573598 0.10164795
## 29.4 0.9934651 0.9926123 0.9971108 0.9985194 0.2185890 0.07545664
## 30 0.9928805 0.9907089 0.9963726 0.9981452 0.1892137 0.05861793
## 30.7 0.9921551 0.9882143 0.9953875 0.9976443 0.1592146 0.04378889
## 31.3 0.9914957 0.9859022 0.9944520 0.9971666 0.1368990 0.03417847
## 32 0.9906823 0.9830877 0.9932803 0.9965651 0.1144526 0.02565550
## group
## Day High - 2 High - 3
## 0 1.000000000 1.0000000
## 0.653 1.000000000 1.0000000
## 1.31 1.000000000 1.0000000
## 1.96 1.000000000 1.0000000
## 2.61 1.000000000 1.0000000
## 3.27 1.000000000 1.0000000
## 3.92 1.000000000 1.0000000
## 4.57 1.000000000 1.0000000
## 5.22 1.000000000 1.0000000
## 5.88 1.000000000 1.0000000
## 6.53 1.000000000 1.0000000
## 7.18 1.000000000 1.0000000
## 7.84 1.000000000 1.0000000
## 8.49 1.000000000 1.0000000
## 9.14 1.000000000 1.0000000
## 9.8 1.000000000 1.0000000
## 10.4 1.000000000 1.0000000
## 11.1 1.000000000 1.0000000
## 11.8 1.000000000 1.0000000
## 12.4 1.000000000 1.0000000
## 13.1 1.000000000 1.0000000
## 13.7 1.000000000 1.0000000
## 14.4 1.000000000 1.0000000
## 15 1.000000000 1.0000000
## 15.7 1.000000000 1.0000000
## 16.3 1.000000000 1.0000000
## 17 1.000000000 1.0000000
## 17.6 1.000000000 1.0000000
## 18.3 0.999999998 1.0000000
## 18.9 0.999999987 1.0000000
## 19.6 0.999999898 0.9999999
## 20.2 0.999999404 0.9999995
## 20.9 0.999995367 0.9999963
## 21.6 0.999964753 0.9999755
## 22.2 0.999804571 0.9998821
## 22.9 0.998630262 0.9993197
## 23.5 0.993256843 0.9972467
## 24.2 0.963270221 0.9882810
## 24.8 0.877693124 0.9671076
## 25.5 0.673350712 0.9168243
## 26.1 0.464647518 0.8515175
## 26.8 0.271226227 0.7586888
## 27.4 0.163985742 0.6736794
## 28.1 0.089673379 0.5759389
## 28.7 0.053310204 0.4971782
## 29.4 0.029136942 0.4136558
## 30 0.017426496 0.3501580
## 30.7 0.009612636 0.2857524
## 31.3 0.005794909 0.2385446
## 32 0.003223455 0.1920658eff.surv3 <- predictorEffect(c("Day"), survmod3a)
surv.effects3<-plot(eff.surv3,
lines=list(multiline=TRUE,
col=c("blue", "blue", "blue", "blue", "red", "red", "red", "red"),
lty=c(1, 2, 4, 3, 1, 2, 4, 3)),
confint=list(style="bands", alpha=0.1),
lwd=4,
axes=list(y=list(lim=c(0, 1.1))),
type="response",
ylab=expression(bold("Prob(Survivorship)")),
legend.position="top",
xlab=expression(bold("Day")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Temperature - Fusion Partners")),
cex=1,
cex.title=1)));surv.effects3
Display mean survivorship within multiple juvenile groups at the end of the experiment.
meanstresssurv2 <- ddply(multiple_final_stress, c("Fusion.Partners", "Treatment"), summarise,
N = length(Survivorship[!is.na(Survivorship)]),
mean = mean(Survivorship, na.rm=TRUE),
sd = sd(Survivorship, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Survivorship, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstresssurv2## Fusion.Partners Treatment N mean sd se max lower
## 1 0 Ambient 9 0.55555556 0.5270463 0.17568209 1 0.02850928
## 2 0 High 26 0.07692308 0.2717465 0.05329387 1 -0.19482341
## 3 1 Ambient 16 0.81250000 0.4031129 0.10077822 1 0.40938711
## 4 1 High 15 0.06666667 0.2581989 0.06666667 1 -0.19153222
## 5 2 Ambient 9 0.77777778 0.4409586 0.14698618 1 0.33681923
## 6 2 High 16 0.00000000 0.0000000 0.00000000 0 0.00000000
## 7 3 Ambient 12 0.83333333 0.3892495 0.11236664 1 0.44408386
## 8 3 High 8 0.00000000 0.0000000 0.00000000 0 0.00000000
## upper
## 1 1.0826018
## 2 0.3486696
## 3 1.2156129
## 4 0.3248656
## 5 1.2187363
## 6 0.0000000
## 7 1.2225828
## 8 0.00000005. STRESS FUSION
Analysis
Analyze with a logistic binomial regression model.
Analyze fusion success over the course of the stress test.
#with full
fusemod2<-glmer(Fusion ~ Day * Treatment * Richness + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family=binomial, link=logit, subset=c(Community=="Multiple"))
summary(fusemod2)
Anova(fusemod2, type="II")
plot(allEffects(fusemod2))
dispersion_glmer(fusemod2) #no overdispersion
plot(residuals(fusemod2)) + abline(h=0, lty=2)
meanstressfuse <- ddply(multiple_final_stress, c("Treatment", "Richness", "Timepoint"), summarise,
N = length(Fusion[!is.na(Fusion)]),
mean = mean(Fusion, na.rm=TRUE),
sd = sd(Fusion, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Fusion, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstressfuse## Treatment Richness Timepoint N mean sd se max
## 1 Ambient 1 End 19 0.7368421 0.4524139 0.10379087 1
## 2 Ambient 2 End 23 0.8260870 0.3875534 0.08081047 1
## 3 Ambient 4 End 4 1.0000000 0.0000000 0.00000000 1
## 4 High 1 End 25 0.6000000 0.5000000 0.10000000 1
## 5 High 2 End 18 0.5555556 0.5113100 0.12051692 1
## 6 High 3 End 14 0.4285714 0.5135526 0.13725270 1
## 7 High 4 End 8 1.0000000 0.0000000 0.00000000 1
## lower upper
## 1 0.28442818 1.189256
## 2 0.43853357 1.213640
## 3 1.00000000 1.000000
## 4 0.10000000 1.100000
## 5 0.04424556 1.066866
## 6 -0.08498116 0.942124
## 7 1.00000000 1.000000Fusion is decreased at the end of the experiment.
6. STRESS GROWTH
Analyze with a poisson regression model.
Analyze growth during the course of the stress test.
Test for effect of temperature:
#with full dataset
#polypmod2<-glmer(Polyps ~ Treatment + Community + Community:Treatment + (1|Parent.Site/Genotype) + (1|Treatment:Tank) + (1|Slide), data=stress, family = poisson, subset=c(Timepoint=="S2D11"))
#summary(polypmod2)
#Anova(polypmod2, type="II")
#qqPlot(residuals(polypmod2))
polypmod2a<-glmer(Polyps^2 ~ Day * Treatment * Community + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family = poisson)
summary(polypmod2a)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: poisson ( log )
## Formula: Polyps^2 ~ Day * Treatment * Community + (1 | Parent.Site/Colony) +
## (1 | Treatment:Tank) + (1 | Slide/Sample)
## Data: stress
##
## AIC BIC logLik deviance df.resid
## 8676.7 8738.9 -4325.4 8650.7 870
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -9.9227 -0.8706 0.0316 0.8371 8.9788
##
## Random effects:
## Groups Name Variance Std.Dev.
## Sample:Slide (Intercept) 0.198862 0.44594
## Slide (Intercept) 0.004889 0.06992
## Colony:Parent.Site (Intercept) 0.191892 0.43805
## Treatment:Tank (Intercept) 0.003544 0.05953
## Parent.Site (Intercept) 0.011778 0.10853
## Number of obs: 883, groups:
## Sample:Slide, 156; Slide, 90; Colony:Parent.Site, 20; Treatment:Tank, 9; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.498522 0.153886 22.735 < 2e-16 ***
## Day 0.010874 0.001113 9.773 < 2e-16 ***
## TreatmentHigh -0.158424 0.139109 -1.139 0.2548
## CommunityMultiple -0.181614 0.124825 -1.455 0.1457
## Day:TreatmentHigh 0.014076 0.002115 6.654 2.85e-11 ***
## Day:CommunityMultiple 0.002843 0.001408 2.019 0.0435 *
## TreatmentHigh:CommunityMultiple 0.265081 0.164704 1.609 0.1075
## Day:TreatmentHigh:CommunityMultiple -0.031383 0.002530 -12.403 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Day TrtmnH CmmntM Dy:TrH Dy:CmM TrH:CM
## Day -0.108
## TreatmntHgh -0.446 0.109
## CmmntyMltpl -0.520 0.133 0.508
## Dy:TrtmntHg 0.049 -0.524 -0.159 -0.060
## Dy:CmmntyMl 0.085 -0.790 -0.086 -0.176 0.415
## TrtmntHg:CM 0.343 -0.093 -0.770 -0.704 0.133 0.127
## Dy:TrtmH:CM -0.041 0.438 0.133 0.089 -0.841 -0.555 -0.171
## convergence code: 0
## Model failed to converge with max|grad| = 0.00180889 (tol = 0.001, component 1)
## Model is nearly unidentifiable: very large eigenvalue
## - Rescale variables?Anova(polypmod2a, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Polyps^2
## Chisq Df Pr(>Chisq)
## Day 302.2053 1 < 2.2e-16 ***
## Treatment 0.8649 1 0.35238
## Community 4.2693 1 0.03881 *
## Day:Treatment 48.8640 1 2.743e-12 ***
## Day:Community 34.2622 1 4.816e-09 ***
## Treatment:Community 0.2662 1 0.60586
## Day:Treatment:Community 153.8404 1 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1qqPlot(residuals(polypmod2a))
## 246 982
## 232 807Plot difference in temperature between single and multiple groups:
stress$group<-paste(stress$Treatment, "-", stress$Community)
polypmod2b<-glmer(Polyps^2 ~ Day * group + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide/Sample), data=stress, family = poisson)
eff.surv5 <- predictorEffect("Day", xlevels=4, polypmod2b)
polyp.effects2<-plot(eff.surv5,
lines=list(multiline=TRUE,
col=c("blue", "blue", "red", "red"),
lty=c(1,2, 1, 2)),
confint=list(style="bands", width=4),
lwd=4,
#axes=list(y=list(lim=c(0, 1))),
type="response",
ylab=expression(bold("(Total Polyp Expansion)"^2)),
legend.position="right",
xlab=expression(bold("Time (Days)")),
main="",
lattice=list(key.args=list(space="top",
border=FALSE,
title=expression(bold("Treatment - Number of Juveniles")),
cex=1,
cex.title=1)));polyp.effects2
Present means of growth by number of juveniles. Sample size was too low due to mortality at the end of the experiment to calculate growth within multiple juvenile groups, so we pool these together for analysis here.
meanstressgrowth <- ddply(stress, c("Day", "Community", "Treatment"), summarise,
N = length(Polyps[!is.na(Polyps)]),
mean = mean(Polyps, na.rm=TRUE),
sd = sd(Polyps, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Polyps, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstressgrowth## Day Community Treatment N mean sd se max lower
## 1 0 Individual Ambient 22 5.727273 2.333643 0.4975343 10 3.393630
## 2 0 Individual High 23 5.000000 1.834022 0.3824200 9 3.165978
## 3 0 Multiple Ambient 42 5.000000 1.847873 0.2851330 9 3.152127
## 4 0 Multiple High 71 5.338028 1.764278 0.2093812 9 3.573750
## 5 3 Individual Ambient 21 6.285714 2.052873 0.4479735 10 4.232842
## 6 3 Individual High 24 5.833333 1.307725 0.2669383 9 4.525608
## 7 3 Multiple Ambient 46 5.891304 1.636096 0.2412293 9 4.255208
## 8 3 Multiple High 63 5.777778 1.754666 0.2210671 10 4.023112
## 9 11 Individual Ambient 20 6.250000 2.022895 0.4523331 10 4.227105
## 10 11 Individual High 19 6.105263 1.370107 0.3143241 9 4.735156
## 11 11 Multiple Ambient 46 6.000000 1.646545 0.2427698 10 4.353455
## 12 11 Multiple High 63 5.904762 1.672494 0.2107144 10 4.232268
## 13 17 Individual Ambient 20 6.500000 2.139848 0.4784844 10 4.360152
## 14 17 Individual High 18 6.166667 1.465285 0.3453709 9 4.701382
## 15 17 Multiple Ambient 45 5.977778 1.815200 0.2705940 10 4.162578
## 16 17 Multiple High 62 5.887097 1.747174 0.2218913 10 4.139923
## 17 23 Individual Ambient 20 7.000000 2.294157 0.5129892 11 4.705843
## 18 23 Individual High 16 6.250000 1.570563 0.3926406 9 4.679437
## 19 23 Multiple Ambient 45 6.577778 1.888829 0.2815701 10 4.688948
## 20 23 Multiple High 57 5.649123 1.894458 0.2509271 10 3.754664
## 21 26 Individual Ambient 20 7.450000 2.480980 0.5547641 11 4.969020
## 22 26 Individual High 4 8.250000 4.349329 2.1746647 12 3.900671
## 23 26 Multiple Ambient 41 6.439024 2.013067 0.3143883 10 4.425957
## 24 26 Multiple High 15 4.266667 1.869556 0.4827172 7 2.397111
## 25 32 Individual Ambient 20 6.500000 2.781518 0.6219663 11 3.718482
## 26 32 Individual High 2 12.000000 0.000000 0.0000000 12 12.000000
## 27 32 Multiple Ambient 36 6.888889 2.213953 0.3689921 10 4.674936
## 28 32 Multiple High 2 1.000000 0.000000 0.0000000 1 1.000000
## upper
## 1 8.060915
## 2 6.834022
## 3 6.847873
## 4 7.102306
## 5 8.338587
## 6 7.141058
## 7 7.527401
## 8 7.532444
## 9 8.272895
## 10 7.475370
## 11 7.646545
## 12 7.577256
## 13 8.639848
## 14 7.631951
## 15 7.792977
## 16 7.634271
## 17 9.294157
## 18 7.820563
## 19 8.466607
## 20 7.543581
## 21 9.930980
## 22 12.599329
## 23 8.452091
## 24 6.136222
## 25 9.281518
## 26 12.000000
## 27 9.102842
## 28 1.000000meanstressgrowth2 <- ddply(stress, c("Day", "Treatment"), summarise,
N = length(Polyps[!is.na(Polyps)]),
mean = mean(Polyps, na.rm=TRUE),
sd = sd(Polyps, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Polyps, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
meanstressgrowth2## Day Treatment N mean sd se max lower upper
## 1 0 Ambient 64 5.250000 2.039296 0.2549121 10 3.210704 7.289296
## 2 0 High 94 5.255319 1.777616 0.1833471 9 3.477703 7.032935
## 3 3 Ambient 67 6.014925 1.770914 0.2163516 10 4.244011 7.785840
## 4 3 High 87 5.793103 1.636345 0.1754346 10 4.156758 7.429449
## 5 11 Ambient 66 6.075758 1.756838 0.2162518 10 4.318919 7.832596
## 6 11 High 82 5.951220 1.601715 0.1768799 10 4.349504 7.552935
## 7 17 Ambient 65 6.138462 1.919285 0.2380580 10 4.219177 8.057746
## 8 17 High 80 5.950000 1.683125 0.1881791 10 4.266875 7.633125
## 9 23 Ambient 65 6.707692 2.013417 0.2497336 11 4.694276 8.721109
## 10 23 High 73 5.780822 1.835200 0.2147940 10 3.945622 7.616022
## 11 26 Ambient 61 6.770492 2.209023 0.2828364 11 4.561469 8.979515
## 12 26 High 19 5.105263 2.941933 0.6749258 12 2.163330 8.047196
## 13 32 Ambient 56 6.750000 2.413974 0.3225809 11 4.336026 9.163974
## 14 32 High 4 6.500000 6.350853 3.1754265 12 0.149147 12.850853Next, analyze the effect of temperature on growth (polyp expansion) at day 23 of 32 during the thermal stress test. As seen in table above, this is the timepoint with sufficient sample size, as later timepoints had high mortality.
First, for the effect of community type (individual vs multiple). No need to square transform as residuals are normal on scale of raw data.
polypmod3a<-glmer(Polyps ~ Treatment * Community + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=d23_stress, family = poisson)
summary(polypmod3a)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: poisson ( log )
## Formula: Polyps ~ Treatment * Community + (1 | Parent.Site/Colony) + (1 |
## Treatment:Tank) + (1 | Slide)
## Data: d23_stress
##
## AIC BIC logLik deviance df.resid
## 607.4 630.8 -295.7 591.4 131
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.0539 -0.4361 0.0267 0.4877 1.6492
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 0.000e+00 0.000e+00
## Colony:Parent.Site (Intercept) 8.515e-03 9.228e-02
## Treatment:Tank (Intercept) 1.351e-10 1.162e-05
## Parent.Site (Intercept) 0.000e+00 0.000e+00
## Number of obs: 139, groups:
## Slide, 81; Colony:Parent.Site, 20; Treatment:Tank, 8; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.94217 0.08946 21.711 <2e-16 ***
## TreatmentHigh -0.11005 0.13286 -0.828 0.407
## CommunityMultiple -0.05117 0.10443 -0.490 0.624
## TreatmentHigh:CommunityMultiple -0.05348 0.15657 -0.342 0.733
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrtmnH CmmntM
## TreatmntHgh -0.614
## CmmntyMltpl -0.794 0.525
## TrtmntHg:CM 0.523 -0.852 -0.665
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(polypmod3a, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Polyps
## Chisq Df Pr(>Chisq)
## Treatment 4.5630 1 0.03267 *
## Community 0.9225 1 0.33683
## Treatment:Community 0.1167 1 0.73266
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1qqPlot(residuals(polypmod3a))
## 380 265
## 53 36There was no difference in growth between individual and multiple juvenile groups, but an effect of temperature was significant. Next, analyze within multiple juvenile groups to look for effect of richness and fusion.
Generate a plot of this relationship.
d23_stress$group<-paste(d23_stress$Treatment, "-", d23_stress$Community)
polypmod3a2<-glmer(Polyps ~ Treatment + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=d23_stress, family = poisson)
eff.stress.polyps1 <- predictorEffect("Treatment", polypmod3a2)
polyps.stress.effects1<-plot(eff.stress.polyps1,
lines=list(multiline=TRUE,
col=c("black"),
lty=1),
confint=list(style="bars", width=4, col="black"),
lwd=4,
axes=list(y=list(lim=c(5,8))),
type="response",
ylab=expression(bold("Polyp Expansion")),
legend.position="right",
xlab=expression(bold("Temperature")),
main=""); polyps.stress.effects1
#lattice=list(key.args=list(space="right",
#border=FALSE,
#title=expression(bold("Genotypic \nRichness")),
#cex=1,
#cex.title=1)));surv.effects1polypmod3b<-glmer(Polyps ~ Treatment * Richness * Fusion.Partners + (1|Parent.Site/Colony) + (1|Treatment:Tank) + (1|Slide), data=d23_stress, family = poisson, subset=c(Community=="Multiple"))
summary(polypmod3b)## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: poisson ( log )
## Formula:
## Polyps ~ Treatment * Richness * Fusion.Partners + (1 | Parent.Site/Colony) +
## (1 | Treatment:Tank) + (1 | Slide)
## Data: d23_stress
## Subset: c(Community == "Multiple")
##
## AIC BIC logLik deviance df.resid
## 458.4 490.0 -217.2 434.4 91
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.08863 -0.48696 0.01669 0.43919 1.42014
##
## Random effects:
## Groups Name Variance Std.Dev.
## Slide (Intercept) 0.000e+00 0.000e+00
## Colony:Parent.Site (Intercept) 1.293e-02 1.137e-01
## Treatment:Tank (Intercept) 0.000e+00 0.000e+00
## Parent.Site (Intercept) 6.070e-10 2.464e-05
## Number of obs: 103, groups:
## Slide, 45; Colony:Parent.Site, 20; Treatment:Tank, 8; Parent.Site, 6
##
## Fixed effects:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.61138 0.35341 4.560 5.13e-06 ***
## TreatmentHigh 0.11566 0.41125 0.281 0.779
## Richness 0.13652 0.21716 0.629 0.530
## Fusion.Partners 0.15557 0.16543 0.940 0.347
## TreatmentHigh:Richness -0.11739 0.23592 -0.498 0.619
## TreatmentHigh:Fusion.Partners -0.19732 0.21296 -0.927 0.354
## Richness:Fusion.Partners -0.06498 0.08519 -0.763 0.446
## TreatmentHigh:Richness:Fusion.Partners 0.06960 0.09835 0.708 0.479
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) TrtmnH Rchnss Fsn.Pr TrtH:R TH:F.P Rc:F.P
## TreatmntHgh -0.840
## Richness -0.943 0.794
## Fusn.Prtnrs -0.842 0.719 0.743
## TrtmntHgh:R 0.853 -0.932 -0.901 -0.675
## TrtmntH:F.P 0.656 -0.817 -0.582 -0.778 0.703
## Rchnss:Fs.P 0.897 -0.759 -0.912 -0.911 0.824 0.711
## TrtmH:R:F.P -0.776 0.871 0.790 0.787 -0.887 -0.904 -0.864
## convergence code: 0
## boundary (singular) fit: see ?isSingularAnova(polypmod3b, type="II")## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: Polyps
## Chisq Df Pr(>Chisq)
## Treatment 3.7823 1 0.0518 .
## Richness 0.0517 1 0.8202
## Fusion.Partners 0.0047 1 0.9454
## Treatment:Richness 0.0793 1 0.7783
## Treatment:Fusion.Partners 0.4499 1 0.5024
## Richness:Fusion.Partners 0.0905 1 0.7635
## Treatment:Richness:Fusion.Partners 0.5008 1 0.4792
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1qqPlot(residuals(polypmod3b))
## 380 862
## 33 74There is no effect of fusion or genotypic richness on growth.
Combine figures for Stress Period:
stress_figs<-plot_grid(surv.effects3b2, surv.effects3, polyps.stress.effects1, labels = c("A", "B", "C"), ncol=3, nrow=1, rel_heights= c(1,1,1), rel_widths = c(1,1,0.8), label_size = 20, label_y=1, align="h");stress_figs
ggsave(filename="Figures/stress_figs.png", plot=stress_figs, dpi=500, width=20, height=6, units="in")7. TEMPERATURE DATA
Load in temperature data from Apex files. Data were recorded every 15 minutes and probes were calibrated to a digital thermometer weekly.
juvtemps<-read.csv("Data/ApexTemps.csv", header=T, na.strings="NA")
#convert Date and Time to date/time format
juvtemps$Date <- as.POSIXct(juvtemps$Date, format="%m/%d/%y %H:%M")
#convert to long format
juvtemps<-juvtemps %>% gather(Tank, Temperature, Tank17:Tank24)
juvtemps$Tank<-as.factor(juvtemps$Tank)
dailymax<-juvtemps[which(juvtemps$TimeDay == "930"),]
dailymax<-juvtemps[which(juvtemps$Phase == "Stress"),]
dailymax$Treatment<-ifelse(dailymax$Tank %in% c("Tank17"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank18"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank19"), "High",
ifelse(dailymax$Tank %in% c("Tank20"), "High",
ifelse(dailymax$Tank %in% c("Tank21"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank22"), "High",
ifelse(dailymax$Tank %in% c("Tank23"), "Ambient",
ifelse(dailymax$Tank %in% c("Tank24"), "High",NA))))))))
maxtemps <- ddply(dailymax, c("Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Temperature, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
juvtemps$Period<-ifelse(juvtemps$Phase %in% c("Stress"), "Stress",
ifelse(juvtemps$Phase %in% c("Stress2"), "Stress",
ifelse(juvtemps$Phase %in% c("Acclimation"), "Stress",
ifelse(juvtemps$Phase %in% c("Hurricane"), "Hurricane",
ifelse(juvtemps$Phase %in% c("Ambient"), "Grow-Out",NA)))))Plot timeseries of temperature data with tank temperature in colors.
#subset temperature measurements
juvtempsHIGH<-subset(juvtemps, Period==c("Stress"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsHIGH$Treatment<-ifelse(juvtempsHIGH$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank19"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank20"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank22"), "High",
ifelse(juvtempsHIGH$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsHIGH$Tank %in% c("Tank24"), "High",NA))))))))
juvtempsHURR<-subset(juvtemps, Period==c("Hurricane"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsHURR$Treatment<-ifelse(juvtempsHURR$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank19"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank20"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank22"), "High",
ifelse(juvtempsHURR$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsHURR$Tank %in% c("Tank24"), "High",NA))))))))
juvtempsGROW<-subset(juvtemps, Period==c("Grow-Out"), c(Date, TimeDay, Phase, Tank, Temperature, Period))
juvtempsGROW$Treatment<-ifelse(juvtempsGROW$Tank %in% c("Tank17"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank18"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank19"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank20"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank21"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank22"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank23"), "Ambient",
ifelse(juvtempsGROW$Tank %in% c("Tank24"), "Ambient",NA))))))))
juvtemps<-rbind(juvtempsHIGH, juvtempsGROW)
juvtemps<-rbind(juvtemps, juvtempsHURR)
TimeSeries1j<-ggplot(juvtemps, aes(x = Date, y = Temperature)) +
geom_line(aes(color = Treatment), size = 1) +
scale_color_manual(values=c("blue", "red")) +
ylab("Temperature (°C)")+
xlab("Date")+
ylim(26,32)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-07-31 10:00:00"))), linetype="solid",
color = "black", size=1)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-08-16 00:00:00"))), linetype="dashed",
color = "black", size=1)+
geom_vline(xintercept = as.POSIXct(as.Date(c("2018-09-20 16:00:00"))), linetype="dotted",
color = "black", size=1)+
theme_classic()+
scale_x_datetime(limits = as.POSIXct(as.Date(c("2018-07-30 09:00:00", "2018-09-20 21:00:00")))) +
theme(legend.position="none")+
theme(text = element_text(size = 18, color="black"))+
theme(axis.text = element_text(size=16, color="black"));TimeSeries1j
Plot by time of day by summarizing the mean temp at each time of day for each tank. Display the mean temperatures of each treatment during the Stress test phase (excluding the pause in treatment due to the hurricane system shutdown).
#subset juv temps for the stress period
juvtempsSTRESS<-juvtemps
juvtempsSTRESS$Period<-ifelse(juvtempsSTRESS$Phase %in% c("Stress"), "Stress",
ifelse(juvtempsSTRESS$Phase %in% c("Stress2"), "Stress",
ifelse(juvtempsSTRESS$Phase %in% c("Acclimation"), "Acclimation",
ifelse(juvtempsSTRESS$Phase %in% c("Hurricane"), "Hurricane",
ifelse(juvtempsSTRESS$Phase %in% c("Ambient"), "Grow-Out",NA)))))
juvtempsSTRESS
juvtempsSTRESS<-juvtempsSTRESS[which(juvtempsSTRESS$Phase == "Stress"),]
QuarterlyJ <- ddply(juvtempsSTRESS, c("TimeDay", "Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N),
max = max(Temperature, na.rm=TRUE),
lower = mean-sd,
upper = mean+sd
)
MeansJtemp <- ddply(juvtemps, c("Period", "Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N)
)
MeansJtemp
MeansJtemp2 <- ddply(juvtemps, c("Treatment"), summarise,
N = length(Temperature[!is.na(Temperature)]),
mean = mean(Temperature, na.rm=TRUE),
max = max(Temperature, na.rm=TRUE),
sd = sd(Temperature, na.rm=TRUE),
se = sd / sqrt(N)
)
MeansJtemp2Generate final figure.
JuvTempsDay<-ggplot(QuarterlyJ, aes(x = TimeDay, y = mean, color=Treatment, fill=Treatment)) +
geom_line(color="black", size = 1) +
#facet_wrap(~Period) +
scale_x_continuous(breaks=c(0, 720, 1440), labels=c("0:00", "12:00", "24:00")) +
scale_color_manual(values = c("blue", "red")) +
scale_fill_manual(values = c("blue", "red")) +
theme_classic() +
geom_ribbon(aes(ymin=QuarterlyJ$lower, ymax=QuarterlyJ$upper), linetype=2, alpha=0.4, color=NA) +
ylab("Temperature (°C)") +
xlab("Time of Day") +
ylim(26,32)+
theme(text = element_text(size = 18, color="black"))+
theme(panel.margin = unit(2, "lines"))+
theme(axis.text = element_text(size=16, color="black")); JuvTempsDay
Combine figures for temperature.
temp_figs<-plot_grid(TimeSeries1j, JuvTempsDay, labels = c("A", "B"), ncol=2, nrow=1, rel_heights= c(1,1), rel_widths = c(1,0.8), label_size = 20, label_y=1, align="h");temp_figs
ggsave(filename="Figures/temp_figs.png", plot=temp_figs, dpi=500, width=14, height=6, units="in")